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Abstract 



We have constructed a complete hydrodynamic theory of nucleation and 
growth in a one-dimensional version of an elastic shear martensitic transfor- 
mation with open boundary conditions where we have accounted for interfa- 
cial energies with strain-gradient contributions. We have studied the critical 
martensitic nuclei for this problem: Interestingly, the bulk critical nuclei are 
twinned structures, although we have determined that the dominant route 
for the formation of martensite is through surface nucleation. We have an- 
alytically solved for the surface nuclei and evaluated exact nucleation rates 
showing the strong preference for surface nucleation. We have also examined 
the growth of martensite: There are two possible martensitic growth fronts, 
viz., dynamical twinning and so-called two-kink solutions. These transforma- 
tion fronts are separated by a dynamical phase transition. We analytically 
derive this phase diagram and determine expressions for the speeds of the 
martensitic growth fronts. 
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1. INTRODUCTION AND MOTIVATION: 



Nature offers a variety of transformations between different crystallographic states. For 
many crystals, the result of such a phase transition involves a dramatic change in the macro- 
scopic shape of the material. Quite often, associated with these macroscopic changes is the 
formation of characteristic domain patterns, sometimes called microstructure. When viewed 
under a transmission electron microscope (TEM), the most commonly observed microstruc- 
tures are those related to twinning, viz., symmetry related variants of the product phase(s) 
oriented in a characteristic pattern. The study of these materials is very desirable because 
many technologically important compounds, such as shape memory alloys and the A15 su- 
perconductors, undergo what are commonly known as martensitic transformations — these 
transformations lead to the formation of the microstructre discussed above 

Martensites are to be distinguished from other materials undergoing structural phase 
transitions by the non-diffusive nature of the transition. To be specific, (i) if the struc- 
tural transformation occurs without a net diffusion of atoms across the crystalline unit 
cell boundary, (ii) the transition is discontinuous (also referred to as first-order), and (iii) 
the transformation involves shear strains, the transformation is said to be martensitic, and 
the resulting material a martensite [0]. These transformations can always be described in 
terms of strains and possibly lattice modulations (or so-called lattice vibrational modes) . 
Here we shall specialize to a subclass of martensitic transitions known as proper ferroelastic 
transitions for which the only relevant variables are the purely elastic strains. 

Martensitic transformations occur by a nucleation and growth processes, and the dynam- 
ical path by which such transformations proceed is the subject of this paper. We explore 
the dynamical evolution of a phenomenological model of a one-dimensional, purely elastic, 
martensitic transformation based on the concept of a strain order parameter. We determine 
both the nucleation and dynamical growth aspects of the transformation for this simple 
model. We also make contact with the formation of microstructure, something that Bales 
and one of us has recently proposed can be associated with the dynamics of martensitic 
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growth. 

The consideration of a lower dimensional model has many advantages: Firstly, many 
of the practical problems {e.g., numerical integration of the equation of motion) associated 
with higher dimensional systems are eliminated. Secondly, lower dimensional systems often 
allow for the possibility of obtaining exact analytical results. We shall indeed find many 
exact analytical solutions, and hopefully these results will be able to be extended into higher 
dimensions in the future. 

Summarizing our main findings: 

I - One may numerically solve for the bulk critical nucleus for such transformations, and 
then obtain an excellent analytical approximation to this structure. The bulk critical nucleus 
is twinned. 

II - For any finite system with open boundary conditions, the critical nucleus may be solved 
for exactly, and we find that it exists at the surface of the system. It's energy is (almost 
exactly) four times less than that of the bulk nucleus (namely, the critical nucleus for a 
system of infinite extent), and it is not twinned. 

III - As shown previously [Q, interesting domain wall motions are associated with the 
growth of martensite when the system is quenched to temperature below the transition tem- 
perature. Here we will show that one can understand the variety of growth fronts that are 
found for these dynamics by analytically deriving interfacial growth speeds for a variety of 
domain walls, and then applying a local stability analysis. 

We believe that apart from the inclusion of thermal fluctuations, something that would 
require simulations based on, e.g., Langevin dynamics, this work represents as complete a 
numerical plus analytical examination of this one-dimensional problem as is possible. 

Our paper is organized as follows. In Section 2 we provide all of the requisite math- 
ematical formalism to understand our one-dimensional model system. Then, in Section 
3 we determine the critical nucleus for both bulk and surface nucleation, with particular 
attention paid to the intimate relation between the bulk and surface solutions. Section 4 
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provides a comprehensive study of the growth dynamics that can occur for initial states 
that are supercritical We also provide analytical work that successfully explains much of 
the observed numerical results, and allows for the exact derivation of various dynamically 
interesting quantities {i.e., the growth speed of a transformation front). We have tried to 
compare, wherever possible, our numerical results to those of relevant experiments. Lastly, 
in Section 5 we summarize our results, and forecast the success of similar phenomenologies 
to the more interesting case of higher dimensional systems. 

2. FORMALISM AND THE EQUATION OF MOTION: 

In this section we will present the formalism necessary to describe a model of a one- 
dimensional (ID) version of a first-order, elastic shear transformation via a Ginzburg- 
Landau (GL) theory. This model was first applied to martensite by Falk and we believe 
that it is the simplest GL potential depending on only one spatial variable that includes 
any of the important characteristics of these transitions. For example, we wish to focus 
on martensitic transitions wherein the high-temperature parent phase has a sufficiently 
high symmetry that the transformation shear strains lead to degenerate martensitic product 
states, something that is a prerequisite for the formation of twinning — we use Falk's GL 
potential since it possesses this feature. 

Falk's mean-field theory of martensitic phase transitions lies in the construction of the 
phenomenological GL free-energy density. To be specific, one considers a displacement field 
u{x,t) and the associated strain e{x,t) = du/dx = d^u. Then, one imagines that the local 
free energy density, of the system can be characterized by a nonlinear function of the 
form: 



In Eq. (|ip. A, B, and C are positive, temperature independent, phenomenological constants 
and 6T = T — T^.. The quantity represents the temperature at which the unstrained 
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(1) 



parent phase (e = 0) becomes unstable. In a linear theory, the undercooling, ST, would be 
related to the linear elastic constant c of the high-temperature parent phase hj c = A 6T. 

The first three terms in the Landau energy density above represent the local response 
of the system to a given strain. The particular symmetry of the local potential, viz. 
/^(e) = e), ensures that both positive and negative shear strains of equal magni- 
tude have the same energy. For an appropriate choice of model parameters, a local free 
energy density similar to the one illustrated in Fig. |l] can be produced. GL densities of this 
form can be attributed to systems for which the (high-temperature) parent phase and (low- 
temperature) product phase(s) are separated by a first-order, elastic shear transformation. 
In our model, temperatures greater than the first-order transition temperature, Ti, lead to 
only one absolute minimum, e = 0, on the energy surface; this strain state corresponds to 
the parent phase. Below Ti, the parent state becomes metastable and two doubly degener- 
ate, stable minima develop at e = ±em (the so-called martensitic strains). These minima 
correspond to symmetry related variants of the product phase and are often refered to as 
martensitic twins. 

We also wish to include so-called non-local elastic forces. For the one-dimensional 
problem that we are considering, the appropriate form of this energy is 



Terms are of this form are responsible for many properties. Firstly, they account for non- 
local forces associated with inhomogeneous strain fields. Secondly, in the context of phase- 
transition theory, such terms represent the domain-wall energy associated with the inhomo- 
geneities of two-phase regions (this is the analogue of the so-called Ginzburg energy found 
in the theory of type-II superconductors). Lastly, these terms break the scale invariance 
implicit in Eq. ([I|) for a bulk system — we shall elaborate on this last feature below. 

The total elastic free energy density is thus /l + Zatl, implying that for a system of length 
L defined by the spatial range ^ < x < ^ the total elastic free energy is 




(2) 



L 




(3) 
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In order for the equation of motion to be specified for this system, one could make 
use of the time-dependent Ginzburg-Landau theory (TDGLT) to construct a first-order in 
time, nonhnear, nonlocal partial differential equation (PDE) for the dynamical evolution of 
the system |0. This approach, however, has been shown to lead to a completely incorrect 
description of the dynamical growth aspects of elastic-shear transformations Instead, 
one must properly account for the hydrodynamic character of the sound waves through the 
inclusion of the kinetic energy density associated with the propagating growth interfaces. 
This approach makes use of the fact that nonlinear elastic models of the type discussed here 
can be treated as nonconvective, hydrodynamic systems 0. 

Propagating disturbances in a solid, e.g., thermal phonons, necessarily involve the dis- 
placement of a finite mass of material and hence contribute to the total mechanical energy 
of the system. Recalling that the displacement field, u{x,t), describes the physical displace- 
ments of the atoms in a solid relative to some chosen undistored system, the kinetic energy, 
T, is given by 

T = \pjl\{d,uf]dx , (4) 

where p is the linear mass density of the undistorted bar. In the spirit of our hydrodynamic 
considerations, we also include a Rayleigh dissipation function, i?, with sound wave viscosity 
7, viz., 

R = \l jl\{d,d,uf]dx , (5) 

that allows for the dissipation of energy of the sound waves It is to be stressed that 
accounting for the damping in this manner ensures that the sound waves are always propa- 
gating at sufficiently long wavelengths; this is the sense in which the system's hydrodynamic 
character is being properly accounted for ||^. 

With these ingredients, the equation of motion plus the boundary conditions may be 
determined precisely as described, e.g., in Ref. [j^ — we stress that we are studying a finite 
system, and thus the boundary conditions are an important part of this problem. The bulk 
equation of motion for the system is 

6 



p{^tu) = {,dlu)[A5T - W^d^uf + hC{d^uf] - D{dtu) + i{dtdlu) , (6) 
and the four boundary conditions are 

A5T{d^u)- B{d^uf + C{d^uf + -i{d^dtu)- D{dlu) = {) ai x = ±f (7) 

and 

{dlu) = at X = ±f . (8) 

By differentiating Eq. once, these equations can also be conveniently restated in terms 
of the shear strain: 

p(92e) = dl[A5Te - Be" + Ce^ - D{dle) + ^{dtc)] , (9) 

with boundary conditions 

{dte) = -{D{dle) - A6Te + Be^ - Ce^) (10) 
7 

and 

(d^e) = 0, at X = ±f (11) 

In certain instances (see below) it will be advantageous to refer to the constitutive equations 
in this form. 

These equations involve a large number (eight) of material parameters. However, we 
can substantially reduce the number of parameters required to model our system through 
scaling analysis, and the result of this analysis for bulk systems has been discussed elsewhere 
PJT7[|. Here we are considering systems of a finite length, and have chosen to rescale our 
dynamical equations through the use of a generalized dimensional analysis for boundary 
valued problems [|10|. The mathematical details of this are summarized in the Appendix A; 
the physics behind the resulting equation of motion and boundary conditions are as follows: 
There are four dynamical units, viz., distance, time, mass and temperature. Rescaling the 
"length" of each of these four dynamical units allows us to eliminate four of the material 
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parameters. Further, we can scale the two parameters depending only on lengths, viz. x and 
the displacement field u, using different scale factors (this is analogous to simply rescaling 
the unit strain), and thus we find that we can reduce our problem down to one involving 
only three material parameters. A, the scaled mass density, 6T, the scaled undercooling, 
and L, the scaled length of the system. The resulting equation of motion and boundary 
conditions are 

A{d^u) = idlu)[ST - 3id,uY + 5id,uY] - (dtu) + (dtdlu) , (12) 
6T{d^u) - (d^u)^ + (d^u)^ + (d^dtu) - (dlu) = 0, at s = ±f , (13) 

(dlu) = 0, at a; = ±f , (14) 

(where in the above system of equations, and from now on, we drop the tildes). 
Further, in terms of the scaled variables the local potential is now given by 

h{x) = l6Te'-\e' + le' (15) 

and thus the martensitic strains are now given by 

e^ST) = ±y . (16) 

Thus, one has the following sequence for the relative stability of the unstrained and marten- 
sitic states 1^: 6T > 1/4, only the unstrained state is locally stable; 3/16 < ST < 1/4, 
the unstrained state is locally stable and the doubly degenerate martensitic states are 
metastable; < 5T < 3/16, the unstrained state is metastable and doubly degenerate 
states are stable; ST <0, only the doubly degenerate martensitic states are locally stable. 

3. THE CRITICAL NUCLEUS: 



In the model discussed above, we have tacitly assumed that the system is coupled to an 
infinite heat bath so that the dynamics can be considered to be isothermal. However, if the 
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system is found in the metastable unstrained state at temperatures such that ST < 3/16, a 
local fluctuation (of the displacement field) can lead to the subsequent formation of a locally 
stable region of martensite which can grow to expel (sometimes) all of the unstrained state. 
That is, the decay of the metastable state {e.g., the unstrained bar) requires an excitation 
called the "critical nucleus" with finite activation energy AE. Any initial displacement 
profile that lies "below" this saddle-point configuration will decay to zero, whereas states 
lying "above" this profile may be able to escape from the basin of attraction and grow to 
the (globally) stable product phase. In this section we shall describe both the novel twinned 
nucleus, which characterizes the bulk nucleus of a spatially infinite system, and the localized 
(in strain) surface states, which correspond to the true critical nuclei (saddle-points) of this 
problem. Then, in the next section we will provide some examples of the fascinating growth 
phenomena that can result as the system approaches steady state. 

We can find the critical nucleus by solving the following nonlinear ODE plus nonlinear 
boundary conditions: 

d^^u - dlui6T - 3{d,uf + 5id,u f) = , (17) 
5T{d^u) - {d^uf + {d,^uf - (dlu) = at x = ±f , (18) 

(dlu) = at X = ±1 . (19) 

The above system of nonlinear equations are derived from considering the static, zero-force 
case of Eqs. (|^,|l^,|l^). The bulk or "saddle-point" critical nucleus that we are interested 
in, to be denoted by Ubcn{x), is the lowest energy, localized configuration satisfying the 
above set of equations — the localization of the displacement field is a consequence of the 
physical constraint that for a bulk nucleus, the boundaries must be unaffected by the spatial 
perturbation. (We stress that this constraint is based on physical considerations, and does 
not arise from a purely mathematical treatment of the formalism presented in the previous 
section.) 
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A. Bulk Solutions: 



In this section we consider a system of infinite extent, thus precluding the possibihty 
that the system nucleates at a boundary. For such a bulk critical nucleus, Eqs. ([T8|) and 



(p!9|) are irrelevant, and the solutions of interest must solve Eq. (12) and satisfy dlu = as 
X ±00 for i = 0, 4. 

The work required to find the critical nucleus is greatly simplified by the observation 
that the equations determining it depend on only one (scaled) material parameter, viz., 6T. 
At this point we choose to refer to the ratio of energies E^/Eb instead of 6T, where Eb is 
the energy barrier separating the unstrained state from the martensitic wells, and Ew is the 
energy difference between the e = and e = ±6^ wells — these two energies are displayed 
in Fig. 1. In terms of 6T 



^(5T)= l-6^T+(l-4^T)v/r^4^T 

Eb^ ' -1 + 65T + (1 - 45T) VI - 4(5T ^ ' 

and from now on we shall refer to this ratio, quite simply because it better distinguishes 
between the differing "shapes" of the local Landau potential Eq. (|l]). Noting that = 
1/6 Ew = Eb, we have chosen to examine Ew/Eb = 0.2, 1.0, and 5.0, corresponding 
to 6T = 0.1823, 1/6, 0.1295. 

We an easily determine one feature of the critical bulk nucleus for this system — symme- 
try arguments plus the physical constraint of a vanishing displacement field at the boundaries 
allow us to specify that the bulk critical nucleus is a localized, symmetric displacement field. 
Namely, since localized symmetric states are lower in energy because there are fewer domain 
walls present, we know that this must be the symmetry of our bulk critical state. 

We have solved for the bulk critical nuclei vs. Ew/Eb numerically. This is a difficult 
problem since apart from the symmetry of the solution, one does not know a priori where in 
function space to begin one's search. In order to minimize this difficulty, we have followed 
an approach based on the full dynamical equations. We assume that the critical nucleus 
should be a configuration approximately of the form 
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u{x) = uo exp(-a;2"/rf2") . (21) 

We consider this functional form because it smoothly interpolates between a classical nucleus 
(n oo) and a Gaussian fluctuation of the displacement field with broad domain walls 
{n = 1) — these two limiting states are shown in Fig. 0. We then specify that d is not 
too narrow, say d = 10 domain wall lengths (see, e.g., Appendix B), and for simplicity take 
n = 1. Thus, for the static displacement field 

u{x) = Uo exp{-x'^/d^) (22) 

at if: = 0, we numerically integrate forward in time under the dynamics of the system to find 
the "critical" value of Uq as a function of E^/Eb, viz., the value of Mq for which perturbations 
above or below will grow or decay respectively. In this way, we have dynamically eliminated 
those initial configurations which could not possibly be candidates for the critical nucleus. 

Having obtained this reduced function space, a numerical solution for the critical nucleus 
becomes tractable: we employ a numerical relaxation method |Tl| that solves the of zero- 
force equation (plus boundary conditions) using as input Eq. ( P^ and the critical value of 
Mo found as described above. These bulk critical nuclei for the above-mentioned ratios of 
Ew / Eb, are shown in Fig. |^. (We stress that these solutions are found for L being sufficiently 
large such that if the length of system is then doubled, our solutions are unchanged — this 
confirms that these are indeed bulk nuclei unaffected by the surfaces.) It is interesting to 
observe that the critical nucleus is a fully twinned strain state. This may be understood 
from Fig. 3b where the strain profile for E^/Eb = 1.0 is displayed. This behaviour is a 
direct consequence of the physical as opposed to mathematical constraint of a vanishing 
displacement field at the boundaries, i.e., a localized, twinned strain state allows for the 
localization of the displacement field. (For example, in any configuration with only one sign 
of strain throughout the entire system, at least one boundary must be displaced.) 

To further check that these numerical solutions behave like saddle-point solutions, we 
have examined their evolution under the full dynamical equations of motion in two ways: 

(i) Ubcn{x) (1 ± e) Ubcn{x) (23) 
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(ii) Ubcnix) -> Ubcn{{l ± e)x) (24) 

where numerically we have used e values as small as 0.001. In the case of (i), it was found 
that the 1 — e configuration was always subcritical, while the 1 + e configuration was always 
supercritical. Further, it was also found that a scaling of the independent variable described 
in (ii) was supercritical (subcritical) if the state was increased (decreased) in size. Thus, 
these critical nuclei have just enough driving force (energetically speaking) and are just big 
enough to drive the system into the martensitic state. 

As the ratio Ew/Eb oo {viz.^ as 5T 0), it becomes increasingly difficult to ob- 
tain converged saddle-point solutions to the zero-force equations. We find that the critical 
nucleus becomes more non-classical, viz., more bell-shaped; eventually, it becomes too nu- 
merically sensitive to relax to a saddle-point solution. However, in the opposite limit of 
Ew/Eb 0, the shape of the critical nucleus approaches that of its classical analogue, viz., 
the step function found from n — > oo in Eq. (pi]) , and for this limit we have found that it 
is possible to produce an excellent approximation to the bulk critical nucleus from purely 
analytical considerations. The details of this analysis are provided in Appendix B, and here 
we summarize our results. 

We consider the zero-force equations now in terms of the strain variable e. Firstly, we 
focus on the equations for an infinite system with the constraint that die = at x = ±oo 
and j = 0, 1, .., 4. One then finds the solutions 



exp[2V(5r(x-Xo)] ^25) 



eoix) = ±2V65T. = = 

\3- 166T + 6 exp[2V6T{x - xq)] + 3 exp[AV6f{x - xq)] 

where xo is some integration constant. Each of these solutions correspond to a strain field 

that is symmetric about x = Xmax where 

1 , /3-165Tn 
Xmax = xo + -^==ln( j (26) 



4VST ^ 3 

Now note that since e{x) is a symmetric function, we know that the boundaries of the system 
must be displaced — put another way, the strain field defined in Eq. is not a localized 
solution. Figure | shows a plot of Eq. (|25|) for Ew/Eb = 1.0. A graphical comparison with 
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the exact localized solutions, viz. those found from the relaxation technique [see Fig. 3b], 
suggest that each of the ±eo(a;) strain configurations correspond to one half of the strain 
strain states that make up the twinned bulk critical nucleus. We have tested this conjecture 
and found that a linear combination of the form 

e{x) = eo{x - xo) - eo{x + Xq) (27) 

(which by the even symmetry of Eq. (^) is an odd function, the same symmetry as our 
numerically determined critical nuclei) does in fact satisfy the bulk zero-force equation 
with a small residual error. This error is a function of xq, a quantity which remains to be 
determined. 

Then, we use this approximation in the following way: if we desire an analytical approx- 
imation to our (exact) numerical result, we can find an xq such that Eq. (pT]) is as close 
as possible to the exact strain profile. Thus, define the residual error to be the integral of 
the square of the difference over the length of the system, and then minimize this quantity 
with respect to Xq. For small Ew/Eb, this procedure leads to superb agreement with our 
(numerically) exact solutions. Figure ^ illustrates the somewhat amazing success of this 
procedure for E\y/Eb = 1.0. At this and lower temperatures we thus find that this function 
yields an excellent fit to the exact solution shown in Fig. 3b. Since we cannot find analytical 
solutions for E\y/Eb ^ 5.0, this procedure is limited to low values of this ratio. 

To understand the success of this approximation, we note that although the nonlinear 
nature of this problem excludes the possibility of any linear combination of zero-force solu- 
tions exactly satisfying Eq. (|1^, such as we have used in Eq. (|27|), the structures that we are 
superimposing on one another are like two solitons. The solitons are particle-like entities, 
and it is usual in soliton theory to ascribe an interaction energy between such particles. 
Usually [0, this interaction decays as an exponential. Thus, in Eq. ( pTj) we have a soliton 
interaction energy that goes as exp(— 2xo). As shown in Fig. 3b, for the ratios of Ew/Eb 
studied here, the value of Xq is quite large, and thus this interaction is very small. This sim- 
ple consideration explains the success of our analytical work in producing an approximate 
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bulk critical nucleus. 



B. A Critical Nucleus at the System's Surface: 

In the previous subsection the bulk critical nuclei were found. We wish to stress that 



in deriving these configurations the boundary conditions given in Eqs. (|T8|,[T9|) were effec- 
tively ignored due to our consideration of a system of infinite extent. Now we include the 
boundaries by focusing on a system that is finite in at least one direction. We will show 
that the bulk critical nuclei are not critical nuclei for these systems — in fact, we shall show 
and analytically derive that the critical nuclei for systems that include at least one free 
surface always possess critical nuclei localized at the surface. Such states are found to have 
a significantly lower energy than the bulk critical nuclei. 

We discovered these solutions quite easily — simply note that localized surface states in 
strain space follow directly from the analytical work on the bulk critical nucleus discussed in 
Appendix B and explicitly stated in Eq. (pSf). Indeed, if the localized, single-humped bulk 
solution eo{x — xo) is centered at one of boundaries, Xmax = the resulting state satisfies 

the boundary conditions and provides a stationary solution to the bulk equation, namely 
Eq. d^). Since such a state is similar to that in shown in Fig. 3b and only involves non-zero 
strains in one fourth of the space that a bulk nucleus exists over, clearly it possesses an 
energy that is (almost exactly) only one quarter of the energy of the critical nuclei solutions 
that we found for the bulk. 

To see that this is indeed a saddle-point solution for a finite system with free boundaries, 
we employ a hydrodynamic approach . As in any hydrodynamic theory with a conserved 



quantity (in this case, momentum), a natural physical interpretation of the equation of 
motion is that of a continuity equation, viz., 

dt{Kdtu) = d^J , (28) 

where J, the one dimensional momentum current density of the system, is given by 
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J = 5T{d,u) - {d,uf + {d,uf + (dtd.u) - (dlu) . (29) 

For the case of static, zero-force solutions, the continuity equation dictates that the (one- 
dimensional) divergence of the momentum current, d^J, must be zero, implying that J is 
a constant. For any solution that is localized in space, such as Eq. (^), it trivially follows 
that in fact J = {e.g., simply consider a region for which the strain is vanishingly small 
and thus the constant J must also vanish). Now note that the nonlinear boundary condition 
in Eq. ( |l8l) is nothing more than the condition that J = at x = ±L/2, and thus this 
boundary condition is satisfied. The second boundary condition, Eq. ([ToD , is satisfied since 
we have placed the maximum of the strain bump at the system's surface. 

The natural question that arises is: which of the above-discussed nuclei is the potent 
nucleus for a large but finite system? To answer this question, recall the common definition 
of the nucleation rate, r, as the ratio of the probability flux j across the saddle point to the 
probability n to be in the metastable well: r = j/n {e.g., see the discussion in Ref. [p^ ). 
We can readily provide the nucleation rate r of a kink at the surface, viz., 

r = S exp{-f3AE/A) , (30) 

where AE is the activation energy for the bulk nucleus. (The prefactor 5* in front of the 
Arrhenius term is a rate constant.) The lower activation energy of a nucleus at the boundary 
compared to that of a bulk nucleus (as mentioned above, this ratio of energies is very nearly 
exactly 1:4) implies that surface nucleation will always be preferred over homogeneous bulk 
nucleation, and thus will certainly be the dominant route to the formation of martensite. 
(Of course, the other common inhomogeneity besides a free surface would be impurities, and 
these would be expected to compete with the surface for the role of "most potent" nucleation 
centre.) 
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4. THE PROPAGATING DYNAMICS OF MARTENSITIC GROWTH FRONTS: 



In this section we will describe and examine the dynamical growth aspects of those 
states which have overcome the nucleation barrier, viz., states that lie "above" the critical 
nucleus, as they approach their steady state profiles. We will thoroughly explain the observed 
dependence of the dynamical evolution of the product phase, viz., propagating martensitic 
growth fronts, on the density A, the undercooling 6T and the length of the system, L — 
this complements the brief outline of results for this phenomenon given in Ref. 0|. In these 
investigations we have numerically obtained the evolution of the system under consideration 
by using a variant method of lines — this provides a high accuracy integration of 
the equation of motion, viz., Eq. @, as well as fully accounting for the time-dependent, 
nonlinear boundary conditions given in Eqs. (^^. 

To be specific, to investigate the dynamical evolution of martensite, we first specify 
the parameters characterizing the system, viz., the scaled temperature, length and mass 
density of the system. (To allow us to focus on a physically relevant range of the scaled 
mass density parameter, we note that recent experimental work on the purely elastic bcc — >■ 
fee transformation of pure Lanthanum [|l^] motivates an initial choice of the A parameter 



of about one |T^. Also, since we are not interested in the dynamics of unstable systems 
we restrict our attention to < 6T < 6Ti.) Then, we choose an initial displacement 
field that is both static and supercritical; e.g., for most of our dynamical studies we used 
the supercritical state defined in Eq. (^) with e ~ 0.01. Finally, we employ our numerical 
integration algorithm [|lj] to follow the temporal evolution of our system. We stress that the 
symmetry and/or profile of the initial state does not qualitatively infiuence the interfacial 
dynamics observed, viz., only the above-mentioned parameters infiuence the qualitative 
aspects of the growth of the product phase. 

For ST just below the first-order transition temperature Ti, the profile in Fig. ^ develops. 
This type of dynamical evolution has been called a "two-kink" growth front — also see Fig. 
1 of Ref. M. Here we utilize the evolution shown in Fig. 6a to develop an approximation 
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for the growth interface found in this temperature regime. This approximation provides us 
with an excellent model from which an accurate theory of the instability of this interface 
can be fully developed. 

Consider the representation of the displacement field shown in Fig. |^. The similarity 
to the growth interface of Fig. 6a for a; ^20 is clear. For simplicity, we have shifted the 
"back end" of the interface to be at x = 0, as well as reflecting the displacement field to be 
positive in the region of interest. We do not possess an analytical expression for this interface 
profile, and thus it is necessary to approximate its functional form. This is accomplished 
by ignoring the "smoothing" imposed by the strain gradients; this leads to the following 
piecewise continuous function for u{x,t) |jl9|: 

, X = 



u{x, t) 



CmX , < X < Vt 

(31) 

+ e2)vt — e2X , vt < X < Vst 



, X > Vst 

In this expression is the martensitic strain, 62 is the strain of the interface, and the speeds 
V and Vs will be described and related below. 

The physical constraint that the interface be coherent at x = Vgt results in the following 
important relation between the speed v and the speed of sound Vg (see below), viz., 

V = - — - — -Vs (32) 

[em + 62) 

This one relationship, plus the analytical evaluation of (already available in Eq. ([16|)), v 
and Vs, will allow for all of the following analysis. 

An illustrative limiting case of Eq. (|32| ) occurs when one considers A ^ 0. This implies 
that the speed of sound Vs becomes infinite, and thus Eq. ( p2D requires that 62 = 0. This 
situation can be understood as the overdamped TDGLT limit of the two-kink solution [Q. 
It corresponds to a kink-type propagating solitary wave of amplitude ±6^ moving to the 
right/left at a speed v. This corresponding analytical solution has already been found by 
Gordon and has the following form: 
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e(a;, t) 



{l + exp[{x-vt)/i]y/^ 



(33) 



where the speed v and £ are uniquely determined by the undercoohng 6T. 

Unfortunately, this solution is entirely unphysical. (This is an important point since 
discussions in the literature have used such a state as a model of martensitic domain wall 
propagation under the name of the Eshelby model.) The finite propagation time of the 
elastic field for the physically relevant scenario of A 7^ immediately suggests that such 
motion must have an unbounded kinetic energy |Q. For example, if one calculates the 
kinetic energy K of Eq. (^) using Eq. (^ for a system of length L one finds that K cc L, 
a result that simply reflects the fact that for this propagating interface the entire system is 
moving! Clearly, this high energy state is never going to be selected as the path through 
which the system approaches its steady state. Thus, the unphysical result contained in 
Eq. (p3D serves to emphasize the importance of including the inertia of the displacement 
field in any accounting of the interfacial dynamics of martensitic growth fronts 0. With 
this in mind, we return to the physically relevant A 7^ problem. 

We wish to describe the interface displayed in Fig. 6a, and to this end we note that the 
complete analytical characterization of this propagating interface requires the evaluation of 
the speed of sound Vg, possibly nonlinear, and the forward kink amplitude 62- Note that in 
the moving interface portion of Fig. 7, Vg characterized the motion in the small strain (62) 
region, and for this reason we evaluate Vs by determining the nonlinear speed of sound. In 
a linear theory, the elastic field propagates at a characteristic velocity determined by the 
undercooling ST and mass density A, viz., 

vl = f (34) 

We require an expression for the nonlinear speed of sound, and to this end we consider the 
"nonlinear corrections" to Eq. ([3^), viz., 




6T + ae^ + he^ 
A 



+ ■ ■ ■ 



(35) 
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where a,b,..., are temperature independent constants and, by symmetry, we have only in- 
cluded the even powers of the strain e. Then, we turn to the local, nondispersive equation 
of motion for the strain field which can be written as 

d^e = ^dl[{5T-e' + e')e]. (36) 

Locally, the term {6T — + e^), approximately corresponds to the nonlinear force constant 
in the parent phase, which we will denote as c. Comparing Eq. (|35|) and Eq. (|36|) we see 
that if a = — 1 and 6 = 1, 

vl = j^l{ST-e' + e') (37) 

represents the 4-th order nonlinear corrections to the speed of sound of an elastic wave of 
amplitude e. In particular, for the nonlinear speed of sound Vs of the strain amplitude 62 
interface, we shall use 



Vs^^j{ST-4 + e^,) . (38) 

The last quantity that we need to determine is the amplitude 62 of the forward kink. 
A physically motivated method for obtaining this quantity relies on energy considerations. 
According to the Lagrangian dynamical formalism P] , the time rate of change of the energy 
E of the system is related to the dissipation function R by 

dtE = -2R . (39) 

Thus, in order for us to derive an expression for 62 we need construct the total energy E of 
the system and the dissipation function R. 

The energy of the system is readily calculated by using our approximation for the dis- 
placement field u{x, t), viz., Eq. (0), and the following expression for the total energy density 
S: 

£ = h{d,uf + ( ^-6Te' - \e' + \e' ) . (40) 
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The total energy is then simply the piecewise spatial integration of the above equation 
over the regions Ji = [0,ft], and I2 = [vt,Vst]. The result of such a calculation yields an 
expression for E, viz., 

E = T^^^ [STe, -lel + U + Utc^ -\el + \el]t, (41) 
(em + 62) 4 3 2 4 b 

and dtE is given by a single time derivative of the above equation. 

Using the same procedure as in [Q, we approximately solve for the dissipation function 

through the use of the analytical expression for the overdamped kink, viz. , Eq. (|33D . At long 

times, the dissipation function can be written as 0, 

R ^ (42) 

where v is the growth speed and £ is a measure of the interfacial width (£ ~ 1 domain wall 
width). The energy balance condition, Eq. (|39|), can now be restated in terms of the only, 
a priori, unknown quantity, namely 62, viz., 



A3Te^--el + -el + -5Te^--el + -el]-^v' = Q, (43) 
(Cm + 62) 4 3 2 4 5 m 

where it should be recalled that Vg and v can both be expressed in terms of 62. We have 

tested this calculation by numerically solving for 62 in Eq. (^) and then using this in 

Eq. (13^ to obtain a value for v. We have found that for the parameters E\y/Eb = 1, 

A = 1, £ = y3/(2e^), V = 0.035, which compares superbly with the value found from the 

numerical integration of the full dynamical equations, viz., Vnum = 0.033. Similar favourable 

comparisons are found for other choices of 6T and A pO| . This excellent agreement between 

the analytical and numerical values for the growth speed v clearly justifies the approximation 

that we have used for u{x,t) in Eq. (|31|). 

One may relate this speed to that in the proper, Ist-order, ferroelastic transition under- 
gone by Lanthanum. One predicts that the two-kink interfacial speed would be of the order 
of 10 % (or less) of the smallest speed of sound, namely about 80 m/s. 

We have found that this interfacial motion does not persist — instead, as the temperature 
of the system is lowered a second type of propagating growth front develops. Figure 8 
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illustrates a typical solution for both the displacement and strain fields. The strain figure 
makes clear the fascinating physics that Bales and one of us observed. This growth front 
has been called "dynamical twin formation" ||^ because the moving interface separating 
product from parent phase leaves behind an alternating structure consisting of both doubly 
degenerate low-temperature variants — these are so-called martensitic twins, and in Fig. 8c 
we display the steady-state profile of a short twinned crystal that results from the dynamics 
shown in Figs. 8a/b. 

Below we shall describe the origin of this interfacial motion. However, here we make one 
straightforward observation contrasting this type of interfacial motion with that of Fig. 
Figure |^ shows a snapshot of the kinetic energy density, K{x, t), for the dynamical twinning 
solution at time t = 360 scaled time units, and illustrates the dramatic localization of the 
kinetic energy density around the interface (x ~ 75) of the transformation front. This is 
clearly different from the entirely delocalized kinetic energy density found throughout the 62 
interface shown in Fig. |T0| (also t = 360 scaled time units) — note the very different scales 
of these two figures. 

Both growth profiles can be understood as arising from the finite propagation time of 
the elastic field found for A 7^ 0. Recall that in our approach sound waves always propagate 
at sufficiently long wavelengths. It then follows that the displacement field far from the 
interface must remain fixed, and thus there is an induced stress in the immediate vicinity of 
the growth front; that is, the parent phase is bent locally at the interface. The magnitude of 
the stress induced in the parent phase depends strongly on the material parameters, A and 
6T, through Eq. (|32|). So, as the growth speed approaches the speed of sound, the 62 strain 
increases and will eventually exceed the spinodal strain determined by d"^ fL^e) /de^ = 0, 
viz.. 



0.1(3- V9-205T) . (44) 

When this occurs the (local) interface becomes unstable and a dynamical phase transition 
arises producing a new type of interfacial motion, viz. that shown in Fig. 8a/b. Summariz- 
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ing, we see that the dynamical twinning shown in Fig. 8a/b is a consequence of the local 
instability of the e2-strain interface of Fig. 6 a/b separating product from parent phase. 
It is worthwhile to re-emphasize that although the t = initial states discussed in this 
paper were twinned nuclei, the dynamical evolutions illustrated in Fig. 6 and Fig. 8 are 
independent of the initial displacement field {e.g., see also Fig. 1 and 2 of where different 
initial conditions were used). 

Based on the above, one can expect that in {A,ST) phase space, a twinning/no-twinning 



phase diagram may be constructed. By setting 62 = Css in Eq. (^3|) and substituting Eq. 
(^) in for V, we can solve for the critical value Ac{ST) [^, viz., 

^1/2 _ ^^m^ssjST — Cg^ + Cgg) I , . 

where it is to be recalled that both and e^g are explicit functions of bT through Eqs. 

viewed in (A,(5T) phase space, this trajectory defines the critical value of A at which 
twinning occurs. This curve is plotted in Fig. ^ along with the results from the numerics. 
The divergence of Eq. (^) when 5T = 0.136 defines the range of undercoolings for which 
twinning will not occur for any value of A. The superb agreement of Eq. (|i5|) with that of 
our numerics clearly justifies our theory of the instability of the growth interface. 



5. SUMMARY AND DISCUSSION: 

We have presented an exhaustive study of the dynamical nucleation and growth in a 
one-dimensional version of a first-order, shear elastic martensitic phase transformation. Our 
approach involves treating nonlinear elastic models of these transitions as nonconvective, 
hydrodynamic systems. We have been able to determine for the first time that the bulk 
critical nuclei of the system are twinned structures. We have obtained approximate analytical 
expressions for these structures, and a comparison of Fig. 3b and Fig. 5 shows that our 
results are in superb agreement with the exact numerical solutions. In addition, we have 
exactly solved for the critical nucleus of a finite, but large system with open boundary 
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conditions, viz., Eqs. ( [T7| , p!8| , [T9D . We have shown that, for this case, the critical nucleus is a 
localized surface state with (almost) exactly 1/4 the energy of the bulk critical nucleus. We 
conclude that surface nucleation is the dominant route for the formation of martensite, at 
least in the absence of impurities. 

A full treatment of the martensitic growth dynamics associated with supercritical nuclei 
has also been presented. A phase diagram connecting the rescaled mass coefficient A and the 
undercooling ST has been derived. We have shown that the transformation fronts propagate 
through the system with a speed given by Eq. (^2|). 

The extension of our work to higher dimensions is in principle relatively straightforward. 
In fact, a large amount of numerical work examining a triangular-to-oblique transition has 
already been completed [0. The formal ingredients are the same as the one-dimensional 
model, except that the numerical analysis becomes much more difficult — the restrictions 
to higher dimensional systems are essentially technical ones. New integration software will 
be needed in order to fully investigate the dynamical evolution of such systems. Analytical 
work will be, obviously, greatly complicated. 

Although extending our one-dimensional phenomenology to completely general 2D and 
3D systems can pose substantial technical problems, we have had much success in looking 
at higher dimensional models with special symmetries. In particular, we have been able to 
use the same methodology as in our ID studies to investigate change of volume first-order 
transitions in a (i-dimensional, nonlinear, nonlocal elastic system {d = 1, 2, oo). In these 
models, the system is reduced down to 1 + 1 dimensions for which the variant method of lines 
and Lagrangian formalism presented in this paper are directly applicable in all (i-dimensions. 
For these systems we find that only surface nuclei exist, and that all dynamical evolutions 
of the system to steady state pass through such nuclei — these results will be presented 
elsewhere. 



23 



ACKNOWLEDGMENTS 



One of us (RJG) wishes to thank the many colleagues who have worked with him on this 
problem over the years: Jim Krumhansl, Baruch Horovitz, Gerhard Barsch, Alan Bishop, 
Jamie Morris, Steve Bales, and Andrew Reid. We also wish to thank Bob Kohn, John Ball 
and Dick James for clarifying comments on competing theories of these transitions. We ap- 
preciate the helpful comments of Ken Vos concerning the solution discussed in Appendix B. 
This work was supported by the NSERC of Canada, and the Advisory Research Committee 
of Queen's University. 



24 



APPENDIX A: DIMENSIONAL ANALYSIS: 



Following the notation of Ref. [T^, let u be the displacement at any position in space, 
and assume that 



u = f{Wi,W2,...,Wio) 



(Al) 



me 

ml 
t'2 ' 



where 

Wi = A, some positive constant of dimension [A] 
W2 = B, some positive constant of dimension [B] 
W3 = C, some positive constant of dimension [C] = [B], 
Wi = D, some positive constant of dimension [D] = 
VF5 = 7, the sound wave viscosity of dimension [7] = 
Wq = p, the linear mass density of dimension [p] = y, 
Wj = 6T, the undercooling of dimension [6T] = 6, 
Ws = L, the system size of dimension [L] = i 
Wq = X, the distance along the bar of dimension [x] = £, 

Wio = t, the elapsed time after an initial strain is applied of dimension [t] = r. 

For this analysis, we have used dynamical units as our fundamental basis, viz., 
Li = i (length) 
L2 = m (mass) 
L^ = T (time) 
L4 = 9 (temperature) 

The corresponding dimension matrix is then simply: 

1 1 1 3 1-10 110 
1 1 1 1 1 10000 
-2 -2 -2 -2 -1 1 
-1 01000 
This matrix has a rank r(B) = 4, whence by the Buckingham Pi-Theorem |TD|, the 
number of measurable dimensionless quantities, under our choice of fundamental units, is 



B 



(A2) 
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k = n — r{B) = 10 — 4 = 6, where n is the total number of independent variables and 
constants appearing in the system. It also follows from the Buckingham Pi-Theorem that 
our system can be re-expressed in dimensionless form where vr = u{Q) is a dimensionless 
dependent variable and Q = {tti, 7r2, Tr^} are dimensionless independent vaiiables and 
dimensionless constants. The exact form of these k dimensionless quantities can be found 
by solving for the null space of the dimension matrix B. The dimensionless displacement 
field u{Q) is found by solving the linear system By = —a, where a is the dimension vector 
of u. Tedious algebra reveals that fl and vr can be written as 



B 

7 



t, 71-2 



■B- 


1/2 


B- 




X, TTs = 




.D. 


3. 



L, TT4 



A 

B 



5T, TTs 



pD 

1^' 



B_ 



(A3) 



and 



TT 



D 



u . 



(A4) 



If we identify tts = A, a rescaled mass coefficient, and rescale tti, ...,7r4 using ttq, we arrive 
at the following dimensionless, rescaled independent variables and constants. 



t 



B' 



t 



X 



TT, 



B 



6 iC 

L, 6T = 7r^7r4 



B 



(A5) 
(A6) 



These scaled quantities, along with Eq. ([A4D, can now be used in Eqs. (|^), (|^ and (^, 
to obtain a complete two-parameter model which depends only on a rescaled undercooling 
and a rescaled mass coefficient (A), viz., 



A{dfu) = {diu)[ 5T - 3{d^uy + 5{diuy ] - (dlu) + (didiu) 



(AT) 



with boundary conditions at x = ± 



2 ' 



(dju) = . 



(A8) 
(A9) 
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APPENDIX B: EXACT SOLUTIONS FOR AN INFINITE DOMAIN: 



In terms of the shear strain, the bulk, zero-force equation that we wish to solve for is 
given by 

die - ( (5Te - + ) = (Bl) 

where the boundary conditions are <9^.e(x) = at x = ±00 and j = 0, 1, 4. 
If we multiply Eq. ( [BTD by dxC and integrate, we obtain 

{d^ef - ( 6Te' - ^e^ + ) = c (B2) 

where c is an integration constant. From the constraint that all derivatives of the strain 
must vanish at ±00, it follows that c = 0. 

Let e = ±l/^/w, so that we are now looking for singular solutions of w. Note that this 
transformation is well defined because e{x) is bounded for all x G (—00, 00). Equation (|B2D , 
with c = 0, then becomes: 

{d^wf - ( 46Tw^ -2w + - ) = (B3) 

3 



Re-writing this as an integral over w and x, we have, 

/ , = / dx (B4) 

Both the left and right hand sides of Eq. ( |B4| ) can be evaluated exactly to give the following 
expression for w and x 

ln(-l + AdTw + ./^V2 - 3w + Q5Tw^) 

— ^^^= = x-Xo (B5) 

where Xq is an integration constant. The inversion of this equation will yield w as a. function 
of X, viz., 



, . 3- 165T + 6 exp[2V5T(a;-xo)] +3 exp[4V5T(a;-xo)] 

w(x) = ^ — B6 

24dTexp[2V6T{x - xo)] 
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Solutions to Eq. are finally obtained by recalling that e = il/y^, whence, 



eo(x — Xq) = ±2\/ 66T^ 



exp[2v5T(x — xq)] 



\3- 166T + 6 exp[2V6T{x - xq)] + 3 exp[AVST{x - a;o)] 



(B7) 



A direct substitution of Eq. (|Bl ) into Eq. (17) confirms that they are indeed zero-force 
solutions. 
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FIGURES 

FIG. 1. A plot of the local elastic free-energy density, fL{e), of Eq. (1) for ^ = S = C = 1 

(corresponding to the scaled units used in our nucleation and growth studies) and 6T = 1/6. The 
barrier height is denoted by Eb and the well depth by Ew- 

FIG. 2. The classical nucleus (long-dash line), corresponding to n — oo in Eq. (21), and a 
Gaussian fluctuation (dot-dashed line) corresponding to n = 1. Again referring to Eq. (21), for 
these figures d = 10 and tto = 1. 

FIG. 3. (a) The bulk critical nucleus found from solutions to Eqs. (17,18,19) for 
Ew/Eb = 0.2, 1/6, and 5 — the tallest curve corresponds to the smallest ratio of Ew/Eb- 
(b) The bulk critical nucleus in strain space for the ratio Ew/Eb = 1.0. 

FIG. 4. A plot of Eq. (25) for Ew/Eb = 1.0. This curve has been shifted so that x = 
corresponds to the Xmax given in Eq. (26). 

FIG. 5. A plot of the analytical approximation to the critical nucleus for Ew/Eb = 1.0 based 
on Eq. (27) using the procedure described in the text for evaluating xq. This strain profile should 
be compared with that shown in Fig. 3b, the (numerically determined) exact critical nucleus for 
this temperature. 

FIG. 6. (a) Evolution of the displacement field for a supercritical nucleus (e = +0.01 in Eq. 
(23)) for Ew/Eb = 1.0. The dashed line is the (static) initial displacement field, and the solid 
lines show how the system progresses at the times t = 0, 180, 220, 260 and 300. (b) Evolution 
of the strain field for the supercritical nucleus studied in (a); is the martensite strain given 
by Eq. (16). The dashed line is the (static) initial strain field, and the solid lines show how the 
system progesses at the times t = 0, 400, 500 and 600. (c) The steady-state t — oo profiles of the 
displacement field (solid line) and the strain field (dashed line). In the context of our model, this 
steady state corresponds to a single twin or "bicrystal" . 
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FIG. 7. An idealized representation of the two-kink growth interface. For the t = state being 
the dashed hne, the sohd hne shows the interface at some later time t. v is the interfacial growth 
speed and Vg is the nonlinear speed of sound, a quantity that we have attempted to evaluate in 
Eq. (38). The exclusion of the nonlocal term, viz., d^u, simplifies the problem to finding only a 

function for u{x,t). 

FIG. 8. (a) Evolution of the displacement field for a supercritical nucleus (e = +0.01 in Eq. 

(23)) for Ew/Eb = 6.5. The dot-dashed line is the (static) initial displacement field, and the 
long-dashed line corresponds to a time of i = 360, and the solid line represents a time oit = 720. 
(b) Evolution of the strain field for the supercritical nucleus studied in (a); is the martensite 
strain given by Eq. (16). (c) The steady-state t ^ co profiles of the displacement field (solid line) 
and the strain field (dashed line) . In the context of our model, this steady state corresponds to a 
polytwinned crystal. 

FIG. 9. The kinetic energy density for the dynamical twinning growth front at t = 360 for 
Ew/Eb = 6.5. For x^75, the kinetic energy density is almost completely localized (see text). 

FIG. 10. The kinetic energy density for the two-kink growth front at t = 360 for E\y/Eb = 1.0. 
For a; ^75, the kinetic energy density is spread out over the 62 growth interface. 

FIG. 11. The renormalized mass density, A, as a function of the undercooling ST. The region 
labeled (/) corresponds to those values of (A, ST) for which twinning will occur. Region {II) will 
only result in the two-kink profile shown in Fig. 6a/b. 5Ti = 0.1875 is the transition temperature 
in scaled units. The solid circles are the numerical data of Ref. [4]. 
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